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1.  Introduction 

The  study  of  many-body  properties  of  ultra-cold  and 
low  density  atomic  samples  is  a  topic  of  considerable 
experimental  and  theoretical  interest.  Substantial 
advances  in  this  new  subfield  of  atomic,  molecular,  and 
optical  physics  have  been  reported  in  the  recent  past,  the 
most  spectacular  one  being  without  doubt  the  observa¬ 
tion  of  Bose-Einstein  condensation  in  weakly  interacting 
alkali  vapors  [1-3].  On  the  theoretical  front,  quantum 
field  theories  of  the  interaction  between  atoms  and  light 
fields  are  being  developed  [4-6],  and  are  finding  appli¬ 
cations  in  the  analysis  of  detection  schemes  for  Bose 
condensates,  nonlinear  atom  optics,  etc. 

While  in  the  experiments  reported  so  far,  Bose- 
Einstein  condensation  has  been  detected  via  the  velocity 
distribution  of  the  atoms  in  the  condensate,  it  is  of  inter¬ 
est  to  examine  alternate  techniques,  such  as  optical 
methods,  which  could  give  clear  and  unambiguous 
signatures  of  the  distinctions  between  a  condensate  and 
an  incoherent  atomic  sample  of  the  same  density  [7-12], 
In  this  paper,  we  use  the  analogy  between  nonlinear 


optics  and  nonlinear  atom  optics  to  propose  and  analyze 
another  measurement  scheme  sensitive  to  the  coherence 
properties  of  the  atomic  sample  under  consideration. 
Based  on  this  analogy,  we  have  predicted  in  earlier  work 
[13]  that  it  should  be  possible  to  achieve  atomic  phase 
conjugation  from  an  atomic  condensate.  Since  phase 
conjugation  is  a  coherent  optical  effect,  which  relies  on 
the  mutual  coherence  between  the  optical  fields  in¬ 
volved,  we  argue  and  demonstrate  numerically  that 
atomic  phase  conjugation  is  rapidly  destroyed  by 
incoherent  processes  such  as  spontaneous  emission. 
Hence,  it  is  not  expected  to  occur  for  an  incoherent 
atomic  sample,  and  is  a  sensitive  probe  of  the  coherence 
properties  and  of  the  lifetime  of  the  condensate. 

Section  2  discusses  our  physical  model,  and  proposes 
an  experimental  scheme  where  a  Bose  condensate  is 
diffracted  from  an  optical  grating  for  atoms  to  produce 
two  counter  propagating  condensate  components,  as 
well  as  a  weak  “probe”  which  is  scattered  off  the 
condensate.  In  a  manner  reminiscent  of  the  optical  case. 
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this  configuration  can  lead  to  matter-wave  phase  conju¬ 
gation  via  four-wave  mixing.  The  basic  approximations 
underlying  our  analysis,  as  well  as  the  nonlinear  master 
equation  which  models  our  system,  are  presented. 

Section  3  reviews  the  coherent  dynamics  of  the  atoms 
and  demonstrates  that  indeed,  phase  conjugation  is 
possible  in  principle.  Because  the  atoms  under  consider¬ 
ation  are  modeled  as  two-level  systems,  though,  the 
theory  of  phase  conjugation  is  somewhat  more  compli¬ 
cated  than  in  the  situation  of  degenerate  four-wave 
mixing  familiar  in  nonlinear  optics,  and  is  more  remi¬ 
niscent  of  the  case  where  the  polarization  of  the  optical 
fields  must  be  taken  into  account.  Spontaneous  emission 
is  included  in  Sec.  4.  In  this  case,  the  nonlinear 
Schrodinger  equation  used  to  describe  the  coherent 
effects  is  replaced  by  a  nonlinear  Hartree-Fock  master 
equation,  which  is  solved  numerically.  Spontaneous 
emission  is  shown  to  drastically  reduce  the  amount  of 
phase  conjugation  available,  thereby  demonstrating  its 
sensitivity  to  the  coherence  of  the  condensate.  In  addi¬ 
tion,  we  compare  these  results  to  the  solution  of  a 
Hartree  nonlinear  Monte  Carlo  wave  functions  numeri¬ 
cal  scheme  that  we  have  recently  proposed,  thereby 
gaining  useful  intuition  about  the  way  a  condensate’s 
description  departs  from  Hartree  as  dissipation  is 
increased.  Finally,  Sec.  5  is  a  summary  and  conclusion. 
Calculational  details  are  relegated  to  appendices. 

In  order  to  avoid  any  confusion,  we  note  now  that  the 
nonlinear  Schrodinger  equation  that  we  use  in  this  paper 
is  not  the  Gross-Pitaevskii  nonlinear  Schrodinger 
equation  familiar  in  the  description  of  Bose  conden¬ 
sates,  but  rather  the  nonlinear  Schrodinger  equation  of 
atom  optics  [4,14],  In  the  former  equation,  the  effective 
nonlinearity  results  from  short-range  interactions 
between  ground  state  atoms,  while  in  the  latter  case  it 
results  from  the  long-range,  near-resonant  dipole -dipole 
interaction  between  ground  and  excited  atoms.  Because 
the  diffraction  grating  involved  in  our  scheme  induces 
real  transitions  between  ground  and  excited  atoms,  this 
latter  interaction  is  dominant,  and  hence  we  neglect  the 
short  range  ground-ground  interactions  here.  After 
having  neglected  this  interaction,  we  do  treat  the  dipole- 
dipole  interaction  as  local,  so  that  our  nonlinear 
Schrodinger  equation  is  itself  local,  just  like  the 
Gross-Pitaevskii  equation.  However,  the  physical 
situation  it  describes  is  quite  different,  and  indeed,  its 
Schrodinger  field  is  vectorial  rather  than  scalar. 

2.  Physical  Model 

In  analogy  to  the  optical  case  of  phase  conjugation  by 
four-wave  mixing,  we  aim  to  realize  a  situation  where  a 
weak  probe  wave  interacts  nonlinearly  with  a  strong 


‘  ‘prnnp’  ’  wave  to  establish  a  grating  from  which  a  coun¬ 
ter  propagating  pump  wave  is  scattered  to  produce  a 
phase-conjugate  wave.  A  possible  way  to  establish  this 
geometry  is  shown  in  Fig.  1,  which  illustrates  how  a 
Bose  condensate  is  diffracted  by  a  matter  waves  grating 
produced  by  a  periodic  evanescent  electromagnetic 
wave. 

,  HI  It, 

I  ;  f  Bose-Einstein 


Fig.  1.  Configuration  of  the  atomic  beams  incident  on  and  diffracted 
from  the  evanescent  light  grating. 

If  the  condensate  impinges  on  the  grating  at  the 
Bragg  (Littrow)  angle,  the  incident  beam  will  be  back- 
diffracted  into  two  components,  one  propagating  in  the 
direction  opposite  to  the  incident  wave  (first-order 
diffraction)  and  the  other,  resulting  from  zeroth-order 
diffraction,  at  some  angle  from  these  “beams.”  We 
already  note  at  this  point  several  important  differences 
between  this  situation  and  the  optical  case:  First,  the 
waves  are  now  de  Broglie  matter  waves,  instead  of  elec¬ 
tromagnetic  fields;  second,  these  waves  have  several 
components,  corresponding  to  all  electronic  states  of 
relevance  for  the  system  at  hand,  so  that  the  fields  under 
consideration  are  vector  fields;  finally,  due  to  two-body 
interactions  between  the  atoms,  their  evolution  is  in¬ 
trinsically  nonlinear,  so  that  no  additional  nonlinear 
“crystal”  is  required. 

In  general,  atoms  impinging  on  the  grating  in  their 
ground  electronic  state  will  be  diffracted  in  a  superposi¬ 
tion  of  ground  and  excited  states,  with  two  major  conse¬ 
quences:  the  spontaneous  decay  of  the  excited  atoms  is 
an  incoherent  process  that  destroys  the  coherence 
between  the  various  interacting  waves,  a  mechanism  that 
is  expected  to  reduce  the  amount  of  phase  conjugation 
achievable.  In  addition,  a  spontaneously  emitted  photon 
can  be  reabsorbed  by  another  atom  in  the  sample.  This 
is  the  physical  mechanism  leading  to  the  near-resonant 
dipole-dipole  interaction,  which  is  the  major  source  of 
nonlinearity  in  the  gedanken  experiment  at  hand. 
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We  have  mentioned  in  the  Introduction  that  the 
dynamics  of  a  Bose  condensate  is  described  by  the 
Gross-Pitaevskii  nonlinear  Schrodinger  equation  [15], 
in  which  the  nonlinearity  results  from  short-range 
two-body  collisions  with  a  scattering  length  of  the  order 
of  the  Bohr  radius.  In  the  presence  of  light  fields, 
however,  the  two-body  potential  is  dominated  at  low 
densities  by  the  dipole-dipole  interaction  between 
ground  and  excited  atoms.  Neglecting  the  short-range 
interactions,  subsequently  approximating  the  dipole- 
dipole  potential  by  a  ‘  ‘contact’  ’  potential,  and  perform¬ 
ing  the  Hartree  approximation  yields  an  effective  single¬ 
particle  nonlinear  Schrodinger  equation  which  forms 
the  basis  of  nonlinear  atom  optics  [4].  Although  this 
equation  has  the  same  structure  as  the  Gross-Pitaevskii 
equation,  it  describes  substantially  different  physics. 

However,  since  spontaneous  emission  cannot  be 
ignored  when  excited  atomic  levels  are  involved  (and  in 
addition  is  the  cause  of  the  dipole-dipole  interaction)  the 
nonlinear  Schrodinger  equation  does  not  give  an  appro¬ 
priate  description  of  the  system.  Rather,  it  must  be  re¬ 
placed  by  a  master  equation  approach.  This  last  equation 
results  from  the  common  procedure  of  eliminating  the 
degrees  of  freedom  of  the  continuum  of  modes  of  the 
electromagnetic  field  in  the  Born-Markoff  approxima¬ 
tion.  Invoking  in  addition  the  Hartree-Fock  factorization 
Ansatz  results  in  the  effective  single  particle  nonlinear 
master  equation  [16,14] 

C>p(“’b)  =-jt  (a  \[H,p]\b)  + (a  \LncP\b) 

f  d{i}[(a  :  l|(V  +  Vc)\2  :  3>(p(2;i)p(3;l) 

+  p(3;  b)p(2;  1)) 

-(p(a;l)p(3;2))  +  p(a;2)p(3;l)) 

<1:2|(V+Ve)|fc:3>].  (1) 

Here  p(a\b)  are  the  matrix  elements  of  the  single 
particle  density  operator  p,  the  “integration  over 
numbers”  such  as  d  { i }  means  a  summation  over  a 
complete  set  of  quantum  numbers  for  the  single  particle 
system,  H  is  the  single  atom  Hamiltonian,  Lncp  is  the 
usual  single  atom  spontaneous  emission  term,  the 
two-body  potential  V  describes  the  dipole-dipole  inter¬ 
action,  and  Vc  is  an  imaginary  potential  which  results 
from  the  elimination  of  the  modes  of  the  thermal  elec¬ 
tromagnetic  field  in  the  Born-Markoff  approximation. 

We  specialize  to  the  case  of  a  system  of  two-level 
atoms  with  lower  electronic  level  |g)  and  upper  elec¬ 
tronic  level  \e).  The  single  atom  Hamiltonian  then  takes 
the  form 


P 

H  =  £ — l-  fj(o0cr+<j_,  (2) 

2m 

where  p  is  the  atomic  center  of  mass  momentum,  m  is 
the  single  particle  mass,  co0  =  ki,c  is  the  atomic  transi¬ 
tion  frequency  and  u±  are  the  usual  pseudospin  raising 
and  lowering  operators.  The  dipole-dipole  potential  V 
takes  the  explicit  form  [17-19] 

fir  i 

V=~2  (AM){y0(*o|ri2|)  -  2  *  (r‘ ^  ' (e  <d) ^ 

x  y2(ko|ri2|)}(er+ ®  cr  +  cr^  ®  cr+)  (3) 

and  the  imaginary  potential  Vc  is 

fir  i 

Vc=i^-(N-l){j0(k0\rn\)-- [l-3(e(r12)-(e(rf»2] 
XJzCkok^l)}]^-®  £T+-  tr+®  o--}  .  (4) 

In  this  expression,  j, ,  y,  are  modified  spherical  Bessel 
functions  of  the  first  and  second  kind  respectively  [20], 
N  notes  the  number  of  atoms  in  the  sample, 
r  =  \d\2kll  2>fr -is  the  spontaneous  emission  rate  and  d 
the  dipole  matrix  element  of  the  atomic  transition.  The 
relative  position  of  the  atoms  is  denoted  r  i2  and  e  (x)  is 
a  unit  vector  in  the  direction  of  the  vector  x.  In  the 
following,  we  approximate  these  potentials  by  a  contact 
potential  of  the  form 

Vo 

v+  V£  =  — -8  (r12)(cr+  ®  a  +  cr  ®  cr+) 

2  U 

iVc  o 

+  — r  §  (r i2)0  ®  (J+  +  a-*  ®  a  ) ,  (5) 

2ko 

an  approximation  adequate  for  ultra-cold  atoms  when 
the  thermal  de  Broglie  wavelength  of  the  atoms  is  much 
larger  than  an  optical  wavelength. 

3.  Coherent  Regime 

Having  established  our  model,  we  briefly  review  mat¬ 
ter  waves  phase  conjugation  in  the  coherent  regime.  In 
this  case,  the  dynamics  of  the  system  is  described  by  the 
Hartree  nonlinear  Schrodinger  equation  [4] 

ih  d^fa)  =  J  d2  (a\H\2)(f>(2) 

+  J  dld2d3  (a,  l|V|2,  3)<£*(l)<^(2)</>(3)  .  (6) 
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Note  that  in  this  equation,  we  have  dropped  the  imagi¬ 
nary  potential  V,  .  While  both  V  and  V ,  find  their 
physical  origin  in  spontaneous  emission,  V  conserves  the 
purity  of  the  state,  while  in  general  Vc  does  not  conserve 
Tr(p2)  [14],  As  such,  it  should  really  be  considered  a 
dissipative  term,  and  consistently,  be  neglected  in  the 
coherent  regime. 

In  complete  analogy  with  conventional  nonlinear 
optics  [21-23]  we  proceed  by  introducing  slowly  vary¬ 
ing  forward  F  (r,  t )  =  ( Fe(r ,  f),  Fg (r ,  t))  and  backward 
B(r,t)  =  (Be(r,t),  Bg{r,t))  propagating  components 
of  the  effective  single  particle  wave  function  </>(r,  t) 

4>(r,  t)  =  F(r ,  +  B  (r,  t)e-i(kz+M\  (7) 

with  \V2  F((r ,  t)\  «  k\V  F  f(r,  t)\  «  k2\F<-(r ,  t)\  and 
|  V2  Bf{r,  t)|  «  k  |  VBe(r,  f)|  «  k2\B  <  (r,  f)|,  where 
£  =  {e,  g }  labels  the  electronic  state  of  the  atoms.  Here, 
we  assume  that  the  main  direction  of  propagation  of  the 
atoms  is  along  the  7-axis,  see  Fig.  1.  When  substituted 
into  the  nonlinear  Schrodinger  equation,  Eq.  (6)  with  the 
explicit  form  of  the  dipole-dipole  potential  V,  this 
decomposition  yields  a  system  of  four  coupled  paraxial 
wave  equations 

(s+2tl)f'=iv;F- 

-iV0(\Fg\2Fe+\Bg\2Fe  +  FgB*Be) 

-iV0(\Fg\2Be+  \Bg\2Be  +  BgF*Fe) 

{i+2k£)F-“vlF- 

-iVQ(\Fe\2Fg+  \Be\2Fg  +  FeB*Bg) 

(s+2tl)s'=iV;B' 

-iV0(\Fe  1 2Bg  +  \Be\2Bs  +  BeF*Fg ) ,  (8) 

where  we  have  introduced  the  dimensionless  wave- 
number  k  — ■>  k/k0,  time  t  — >  wre cf,  wavefunction 
<t>(y,z,t)^><t>(y,Z,t)lk0  and  position  r  kQr 
variables,  as  well  as  a  dimensionless  potential 
Vo  — »  yo/wrec .  Equations  (8)  resemble  the  usual  model 
equations  for  beam  propagation  in  nonlinear  optics  [21], 


except  that  due  to  the  vectorial  character  of  the  slowly 
varying  amplitudes  F(r,t)  and  B(r,  t)  and  the 
exchange  of  excitation  involved  in  the  dipole-dipole 
interaction,  the  nonlinear  term  only  couples  partial 
waves  corresponding  to  different  electronic  states.  As 
such,  these  equations  are  reminiscent  of  those  involving 
polarization  coupling  in  nonlinear  optics,  although  in 
the  latter  case,  terms  involving  self-phase  modulation 
would  also  appear.  Still,  the  analogy  with  nonlinear 
optics  is  sufficient  in  that  we  can  readily  adapt  the  beam 
propagation  method  common  in  optics  to  the  problem  at 
hand. 

For  numerical  purposes,  we  restrict  our  attention  to 
two  dimensions  and  to  the  portion  of  the  condensate 
confined  by  the  planes  z  =  0  and  z  =  L  and  assume  that 
the  boundary  conditions  are  provided  by  the  waves 
incident  on  these  planes  from  the  outside  of  that 
volume.1  If  we  assume  that  the  Bose  condensate  consists 
of  ground  state  atoms  sufficiently  far  from  the  grating, 
we  have 

Ft(y,  0;  t)  -  ^0(0,  t)  +  3U0,  t)e,Klf 

Fg(y ,  0;  t)  =-  ^s0(0,  t)  +  &gs(0,  t)eIKy 

Be(y,L ;  0  =  0 

Bg(y,L;t)=dBt0(L,t)  .  (9) 

Here,  we  have  introduced  a  spatial  Fourier  decomposi¬ 
tion  of  the  transverse  variation  of  the  amplitudes  to 
account  for  diffraction  by  the  grating,  which  is  assumed 
to  leave  the  atoms  in  a  superposition  of  their  excited  and 
ground  electronic  states.  Hence  (z,  0 ,  S'go (z,  0  and 
S8so  (z,  t)  are  the  amplitudes  of  the  forward  and  back¬ 
ward  propagating  pump  beams  along  the  7-axis,  while 
3Fes(z,  0  and  S'gaiz,  0  are  the  amplitudes  of  the  forward 
propagating  probe  waves  propagating  at  a  small  angle  to 
the  7-axis  such  that  k«  k. 

For  pump  waves  strong  enough  that  the  undepleted 
pump  approximation  holds,  i.e.,  9m(z,  0  —  const, 
3ftfo(z,  0  ~  const,  this  problem  can  be  solved  analyti¬ 
cally  [21,22],  It  was  shown  in  Ref.  [13]  that  the  intrinsic 
nonlinearity  of  the  system  results  in  the  ground  state 
probe  being  associated  to  an  excited  state  conjugate 
wave  propagating  with  transverse  wave  number  -k  ,  and 
similarly,  the  excited  probe  leads  to  a  ground  state 
conjugate. 


The  .v -dimension  is  eliminated  by  expressing  the  wave  function 
<f>e(r,t)  as  <b(x)<b/(y .  z.,  I )  where  <l> (x )  is  taken  to  be  a  solution 
of  a  linear  Schrodinger  equation.  As  a  result  Vo  is  replaced  by 
Vo  f  dx\<f>(x)l4/ko. 
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We  performed  numerical  simulations  of  the  coupled 
paraxial  wave  equations,  Eq.  (8)  for  various  boundary 
conditions  involving  different  fractions  17  of  excited 
atoms  in  the  diffracted  wave  at  7  =  0.  As  a  numerical 
check,  we  verified  that  the  atomic  flux  through  the  z  =  0 
and  z  =  L  planes  is  properly  conserved.  Fig.  2  shows  an 
example  of  the  build-up  of  the  ground  state  phase  conju¬ 
gate  wave  of  amplitude  Sftg(-  k,  0,  f),  and  the  associated 
depletion  of  the  pump  amplitude  3ft g  (0,  0,  t). 

P 


Fig.  2.  Time-dependence  of  the  backward  propagating  atomic 
beams|  58g(-  k,  0.  f ) | 2  —  dashed  line,  and  1 58g(0,  0,  f)|2  —  solid  line, 
for  K  =  k0,  k=mo,L  =  40/ko,  Vo  =  0.4o>rec,  VrO  =  0.  T  =  0, 
=  0.5,  =  0.5.  3?g„  =  0,  9SS  =  0  and  S8g0  =  0.7. 

The  dipole-dipole  potential  V  responsible  for  the  non¬ 
linearity  is  switched  on  at  time  t  =  0  and  the  atoms  are 
assumed  to  be  left  in  their  excited  state  after  being 
diffracted  by  the  grating.  In  that  case,  the  conjugate 
beam  builds  up  only  in  the  electronic  ground  state  as 
expected  from  the  analytical  solution  in  the  undepleted 
regime  [13].  The  partial  Schrodinger  waves  first 
undergo  a  transient,  with  discontinuities  in  their  slope 
after  the  beams  cross  the  condensate  an  integer  number 
of  times.2  The  kicks  become  less  and  less  pronounced  as 
a  new  steady  state  corresponding  to  V  ¥=  0  is  established. 

Figure  3  shows  how  the  fraction  of  excited  state  atoms 
in  the  reflected  beam  influences  the  intensity  of  the 
conjugate  beams.  As  follows  both  from  analytical  con¬ 
siderations  [13]  and  from  our  numerical  simulations,  no 
conjugate  beams  are  produced  if  the  atoms  exit  the 
grating  in  their  ground  state.  This  is  simply  because  in 
that  case,  there  is  no  spontaneous  emission,  and  no 
dipole-dipole  interaction  between  the  atoms.  Note  the 
asymmetry  between  the  excited  and  ground  state  atomic 
conjugation,  which  results  from  the  asymmetry  in  the 


2  This  is  an  artifact  of  our  numerical  algorithm,  which  starts  propagat¬ 
ing  the  states  from  the  two  boundaries  at  z  =  0  and  z-L  after  V  is 
switched  on. 


pump  beams,  the  condensate  being  assumed  to  consist 
of  ground  state  atoms.  Of  course,  the  results  for  large 
fractions  of  excited  atoms  should  be  regarded  with  cau¬ 
tion,  due  to  the  neglect  of  the  incoherent  effects  of 
spontaneous  emission  in  these  results.  This  is  the  object 
of  the  next  section. 


Fig.  3.  Strength  of  the  conjugate  beams  1 (-k  ,  0,  r)  | 2  —  dashed 
line,  and  |  SSg  (-K,  0.  t)\2  —  solid  line,  in  steady  state,  as  a  function  of 
the  fraction  17  of  excited  atoms  in  the  reflected  beam.  K  =  ko,k=  10Ao, 
L  =  40/ko,  To  =  0.4ft) „c,  Vc0  =  0  and  T=  0. 

4.  Spontaneous  Emission 

4.1  Master  equation  approach 

In  order  to  account  for  the  effects  of  spontaneous 
emission,  the  Hartree  nonlinear  Schrodinger  equation, 
Eq.  (6)  must  be  replaced  by  the  nonlinear  master 
equation,  Eq.  (1),  as  already  mentioned.  In  analogy  with 
the  coherent  case,  we  now  expand  the  effective  single 
particle  density  matrix  element  pec(r,  r',  t )  as 

pw{r,r',t)  =  pfrir.r',  t )ett(z'z,)  +  p$(r,r’, 

+  p&  (r,  r',  t)e‘Hz+z,)  +  p%  (r,  r',  t)&ik(z+z,) .  (10) 

This  decomposition  yields  16  coupled  paraxial-like 
partial  differential  equations  of  the  spatial  variables 
( y,z ),  (y',z')  and  1.  Their  explicit  form  is  given  in 
Appendix  A.  These  equations,  being  first-order  in  z  and 
z' ,  require  boundary  conditions  on  the  four  sides  of  the 
square  z,  z’  e  [0,  L],  However,  the  boundary  conditions 
Eq.  (9)  on  the  wave  functions  in  the  coherent  case  yield 
boundary  conditions  for  the  density  matrix  at  the  cor¬ 
ners  of  this  square  only.  The  additional  boundary  equa¬ 
tions  required  for  the  sides  of  the  square  can  be  derived 
from  the  Hartree  nonlinear  Schrodinger  equation  with 
the  complex  potential  V  +  V, ,  as  shown  in  Appendix  A. 
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We  have  numerically  solved  the  Hartree-Fock  nonlin¬ 
ear  master  equation  in  two  dimensions  for  various  spon¬ 
taneous  decay  rates  and  fractions  77  of  excited  atoms  at 
the  exit  of  the  beam  splitter.  As  we  are  interested  in 
particular  transverse  directions  of  the  wave  propagation 
we  introduce  a  spatial  Fourier  transform  of  the  atomic 
density  operator  in  the  transverse  dimension 

p“/  (ky,  z;  k'y,  z',  t)  = 

jdkydk'yP$(y,z;y',z',  t)eik>y  eik'yy,  (11) 

where  a  and  /3  stand  for  the  forward  and  backward 
propagating  wave.3 

Figure  4  shows  the  probabilities  p  “  (-k,  0;  -#c,  0,  t ) 
and  p  "  (0,  0;  0,  0,  t)  as  a  function  of  time,  for  the  case 
where  the  atoms  exit  the  diffraction  grating  in  their 
excited  state.  As  such,  it  is  the  same  as  Fig.  2,  except  that 
now  the  spontaneous  emission  rate  has  been  set  equal  to 
r  =  0.5  and  the  imaginary  potential  is  V ,  =0.1,  in  recoil 
units.  Clearly,  the  amount  of  phase  conjugation  is  sub¬ 
stantially  reduced  as  compared  to  the  coherent  situation. 
This  is  further  illustrated  in  Fig.  5,  which  shows  the 
steady  state  phase  conjugate  probability  as  a  function  of 
the  spontaneous  decay  rate. 

Already  for  a  decay  rate  equal  to  the  recoil  frequency, 
phase  conjugation  has  all  but  disappeared.4 


P 


Fig.  4.  Master  equation  simulations  in  the  incoherent  regime.  Time- 
dependence  of  the  backward  propagating  atomic  beams 
Pgg( 0,  0;  0,  0,  t)  —  solid  line,  p  gg  (-/c,  0;  —k,  0;  t) — dashed  line,  for 
k  =  k(h  k  =  10&o,  T  =  40/&o,  V0  -  0.2o)rec,  Vc0  =  0Aa)rec,  r  =  0.5  corec, 
&spont  =  k0,  &eo  =  0.5,  =  0.5,  &go  -  0,  SFgs  =  0  and  ^g0=0.7. 


4  The  parameters  of  our  simulations  were  largely  dictated  by  com¬ 
puter  limitations.  It  is  expected  that  phase  conjugation  will  survive 
larger  spontaneous  emission  rates  if  larger  dipole-dipole  potentials  and 
longer  propagation  distances  are  considered. 

3 

As  it  should  be,  our  numerical  results  agree  with  the  coherent  ones 
for  r = Vc  =  0. 


Fig.  5.  Influence  of  /  on  the  strength  of  the  conjugate  beam 
(-k,  0; -k.  0;  t)  in  steady  state  K  =  k0 ,  k=l0ko,  L  =  40/ko , 
Vo  0.2oj,„  .  V cO  0.  I  OJj/.f  ,  /cspont  /((,,  .'A Q  0.5,  -iA  =  0.5,  0, 

9gs  =  0  and  5Fso  =  0.7. 


4.2  Hartree  Monte  Carlo  Simulations 

The  numerical  solution  of  the  nonlinear  master  equa¬ 
tion  leads  to  substantial  computer  memory  require¬ 
ments.  In  a  recent  paper  [24],  we  have  proposed  an 
alternative  approach  to  its  brute  force  solution,  based  on 
the  average  over  ‘  ‘quantum  trajectories’  ’  resulting  from 
the  solution  of  a  Hartree  nonlinear  Schrodinger  equa¬ 
tion,  interrupted  by  quantum  jumps.  The  validity  of  this 
approach,  which  is  inspired  by  the  now  well-established 
Monte  Carlo  wave  functions  technique  [25-27],  can 
however  only  be  proven  provided  that  the  Hartree  ap¬ 
proximation  is  exact,  which  holds  strictly  only  for  con¬ 
densates.  In  addition,  the  two-body  contribution  to  dis¬ 
sipation  must  be  of  the  form  of  a  complex  potential, 
which  is  the  case  in  the  present  situation.  A  direct  com¬ 
parison  of  the  Hartree  Monte  Carlo  and  master  equation 
results  should  therefore  give  us  not  only  an  indicator  of 
the  practical  validity  of  the  former  method,  but  also  a 
sense  of  how  spontaneous  emission  destroys  the  conden¬ 
sate. 

In  order  to  implement  this  approach,  we  generalize 
the  Hartree  nonlinear  Schrodinger  equation,  Eq.  (6)  to 
include  the  imaginary  potential  Vc,  so  that  it  becomes 

ifi  6  ^({a>  =  jd2(a\H\2)  (f)(2) 

+Jdld2d3  (a,  l\V+  Vc\2,  3)  <f>*(  1)  0  (2)  (/>  (3) . 

(12) 

Expressing  furthermore  the  linear  Liouvillian  L„c  in 
Eq.  (1)  as  the  Lindblad  form 
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Lnc  p  (t)  =  -  i  2  {cm  cm  p{t)  +  p  (t)  Cl  C,„) 

+  2  Cmpit)Cl  (13) 

m 

where  the  Cm’s  are  system  operators  appearing  in  its 
coupling  to  the  continuum  of  electro-magnetic  field 
modes,  the  Hartree  Monte  Carlo  algorithm  proceeds  as 
follows  [24],  We  first  formally  reexpress  the  nonlinear 
Schrodinger  equation,  Eq.  (12)  as 


In  case  no  jump  occurs,  the  state  vector  evolves  under 
the  influence  of  the  nonlinear  and  nonhermitian  effec¬ 
tive  Hamiltonian  He ff,  so  that 

m,  iffeff  [|<KO>]|<KO>  „ 

(f  +  Sr))  —  \<f>(t))~ - t -  8 1 

n 

(18) 


and  the  new  wave  function  is 


ih  I  (MO)  =  H„,[  |<«1)>]  1 0(0).  (14)  \4>(t  +  8f))  =  |</>U)  (r  +  8f))/VN8/? .  (19) 


where  the  nonlinear  “Hamiltonian”  H„i  [\<f>  (t))] 
is  a  functional  of  the  effective  single-particle  state 
vector  1 0(f)),  see  Eq.  (12).  We  then  decompose  the 
Liouvillian  Lnc  into  a  contribution  //,,  added  to  //„,  to 
yield  the  effective  nonhermitian  “Hamiltonian” 

H cff  [|  <^>)]  =  H„i  [ | cf>)]  +  Hti 

=  H„,[ |<A>]-f  E  Cl  Cm  (15) 

m 

and  a  linear  “fill-up  term,”  see  Eq.  (13).  As  in  the 
usual  case,  Heff  is  responsible  for  the  Schrodinger-like 
evolution  of  the  system,  while  the  “fill-up  term”  causes 
randomly  distributed  quantum  jumps. 

The  probability  8 p  for  one  of  these  jumps  to  occur  is 
given  for  the  small  time  interval  8 1  by5 


/8f  ( (ft(0l//d-HJl<H0) 
A  (<£(0|</>(0) 


<«K0l2  clcMit)) 

- — i - =  2  8 p, 

m) \m)  ™ 


(16) 


We  performed  Monte  Carlo  simulations  within  the 
slowly  varying  envelope  approximation  introduced  in 
Sec.  2  and  for  the  same  beam  configuration  as  in  the 
proceeding  example  (see  Appendix  B).  As  pointed  out  in 
Ref.  [24],  the  accuracy  of  the  Monte  Carlo  simulations 
as  compared  to  the  Hartree-Fockmaster  equation 
description  is  expected  to  depend  on  the  amount  of 
dissipation  in  the  problem.  Figure  6  shows  the  relative 
difference  between  the  two  methods  as  a  function  of  the 
spontaneous  emission  rate  F .  These  results  confirm  that 
the  Hartree  Monte  Carlo  simulations  work  surprisingly 
well,  although  the  discrepancy  between  the  two 
methods  increases  with  F,  as  expected. 


A  jump  is  said  to  occur  if  a  quasi-random  number  e, 
uniformly  distributed  between  0  and  1,  is  less  than  8 p. 
The  post-jump  wave  function  is  chosen  amongst  the 
various  possible  final  states  Cm\(f>(t))  according  to  the 
probability  law  8p,„/8p,  with 

8f)>  =  Cm\<l>(t))/(bpm/?>t)m  .  (17) 


5  Note  that  we  keep  track  of  the  norm  of  the  state  vector  explicitly, 
since  only  a  partial  volume  of  the  system  is  considered  explicitly,  and 
thus  the  flux,  rather  than  the  norm,  is  conserved  in  our  numerical 
solution. 


Fig.  6.  Comparison  of  the  Monte  Carlo  and  master  equation  simula¬ 
tions.  Influence  of  Con  A,  the  standard  deviation  between  the  steady 
state  results  of  the  Monte  Carlo  and  master  equation  predictions, 
normalized  by  the  Master  equation  result. 

We  attribute  this  surprisingly  good  agreement  to  the 
fact  that  the  type  of  dissipation  that  we  are  considering 
here  is  however  quite  special,  in  that  the  only  effect  of 
the  nonlinear  Liouvillian  Z„,  is  in  the  form  of  an  imagi¬ 
nary  two-body  potential  Vc .  It  corresponds  to  a  modifi¬ 
cation  of  the  spontaneous  decay  rate  resulting  from  the 
interference  of  the  spontaneous  emission  probability 
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amplitudes  due  to  the  undistinguishability  of  the  two 
particles.  In  other  words,  it  describes  a  collective  decay 
mechanism  (albeit  as  a  result  of  the  Born-Markov 
approximation,  only  two-atom  effects  are  retained).  Due 
to  this  collective  nature,  it  is  reasonable  to  expect  that 
this  imaginary  potential  might  not  do  too  much  violence 
to  the  Hartree  ansatz. 

The  only  other  incoherent  contribution  to  the  nonlin¬ 
ear  master  equationis  a  single-particle  Liouvillian  L, . 
But  due  to  the  indistinguishability  of  the  particles  in  the 
sample,  it  describes  a  decay  mechanism  that  occurs  for 
all  N  atoms  simultaneously.  As  a  result,  it  leaves  all 
atoms  in  the  same  state,  albeit  now  described  by  a 
product  of  N  identical  single-particle  density  operators. 
It  appears  then  that  for  the  very  specific  form  of  decay 
we  are  considering,  the  Hartree  ansatz  might  not  be  as 
bad  as  otherwise  expected,  a  conjecture  largely  con¬ 
firmed  by  the  result  of  Fig.  6. 

5.  Summary  and  Conclusions 


6.  Appendix  A.  Master  Equation 
in  the  Slowly  Varying  Envelope 
Approximation 

In  analogy  to  the  wave  function  description  in  terms 
of  forward  and  backward  propagating  beams  each  of  the 
four  density  matrix  elements  pee ,  pgg ,  peg ,  pge  can  be 
decomposed  into  four  parts 

pec  (r,  r',  t )  =  pfc  ( r ,  r\  /)e'fc(z“z'J  +  p$  (r,  r',  r)e_'*(z-z') 
+  p{bc  (r,  r',  t)e‘Kz+z,)  +  p$  (r,  r',  t)eik(z+z,) ,  (20) 

where 

p&  =  Fe(r,t)B*(r',t).  (21) 

The  16  coupled  paraxial-like  partial  differential  equa¬ 
tions  obtained  from  the  nonlinear  master  equation, 
Eq.  (1)  in  the  slowly  varying  envelope  approximation  are 
of  the  form 


In  this  paper  we  have  discussed  the  possibility  of 
observing  atomic  phase  conjugation  off  Bose  conden¬ 
sates,  and  to  use  it  as  a  diagnostic  tool  to  determine  the 
spatial  coherence  properties  and  lifetime  of  that  conden¬ 
sate.  Since  phase  conjugation  results  from  the  scattering 
of  a  partial  matter  wave  off  the  spatial  grating  produced 
by  two  other  waves,  it  offers  a  very  natural  way  to 
directly  measure  these  coherence  properties,  and  as 
such  provides  an  attractive  alternative  to  the  optical 
methods  proposed  in  the  past. 

As  expected,  incoherent  mechanisms  such  as  sponta¬ 
neous  emission  rapidly  destroy  phase  conjugation, 
thereby  demonstrating  that  the  proposed  technique 
should  easily  distinguish  a  condensate  from  an  incoher¬ 
ent  sample  of  the  same  density.  We  have  also  compared 
Hartree-Fock  master  equation  results  to  those  based  on 
Hartree  Monte  Carlo  wave  function  simulations,  thereby 
obtaining  a  measure  of  the  destruction  of  the  condensate 
as  a  function  of  the  rate  of  spontaneous  emission. 

We  conclude  by  noting  that  in  replacing  the  nonlocal 
dipole-dipole  potential  by  a  local  interaction,  the  grat¬ 
ings  established  in  the  condensate  are  probed  only 
locally.  It  will  be  interesting  to  investigate  if  more  infor¬ 
mation  can  be  gained  by  treating  V  as  a  nonlocal  poten¬ 
tial.  In  particular,  one  can  expect  that  effects  such  as 
spatial  diffusion  in  the  sample  should  be  accessible  by 
such  nonlocal  sampling. 


(_d_  _d_ 
\dz  +  dz' 


pll(r,r';  t)  = 
i  (Vi  -V±2.)  p//(r,  r'\  t) 


~  iV*  [(pi (r,  r ;  t )  +  p“  (r,  r ;  t))  pfg  (r,  r'\  t)  + 
P%(r,r;t)pg(r,r’;t)] 

~  iVm  Kp/* t)  +  Pgg  0)  pll  (r,  r'\  t )  + 

P$(r,r;t)  p%(r,r';t)] 

+  iVm[(pg(r',r’;  t )  +  p“  (r',r't))  plfg(r,r'\  t)  + 
Pgl(r',r'\t)  pi(r,r'\t)] 

+  iv,0[(plfg  <T,  r'\  t )  +  p“  (r\  r'\  t))  p// (r,  r’;  t)  + 
Pgg(r',r';t)  pS  (r,r';f)] 

+  (er'\Lncpll(t)\er),  (22) 


where  V,0  =  Vo  +  Vo ■  Similar  equations  can  be  obtained 
for  all  other  components. 

These  first-order  differential  equations  in  z  and  z' 
require  boundary  conditions  on  the  four  sides  of  the 
square  z,  z'  =  [0,  /,].  However,  the  boundary  conditions 
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for  the  wave  functions  in  the  coherent  case  Eq.  (9)  yield 
boundary  conditions  for  the  density  matrix  only  at  the 
corners  of  the  square.  Ignoring  spontaneous  emission  at 
the  boundaries,  we  can  express  the  density  matrix 
elements  as  p£  (y,  z,  y',  0,  t)  =  Fe  (y,z,  t)F*(y',  0,  t ), 
etc.  In  this  expression,  Fe(y,z,t)  undergoes  a 
Schrodinger  evolution,  while  F*(y',  0,  t)  is  a  boundary 
value.  Thus,  the  dynamics  of  the  system  on  the  sides 
of  the  square  is  obtained  by  considering  a  boundary 
equation  for  the  master  equation,  which  involves  only 
the  evolution  corresponding  to  the  change  of  “half’  of 
the  density  matrix  element. 

This  boundary  equation  can  be  derived  from  the 
Hartree  nonlinear  Schrodinger  equation  with  the  poten¬ 
tial  V  +  Vc.  In  particular,  we  find  for  pi  ( r ,  r'b0und;  0 


resembles  atomic  recoil  due  to  the  spontaneous  emission 
event  [26].  The  effective  Hamiltonian  He ff  [|  (/>  )]  is  then 

ifiF 

Helt[\cl>)]  =  Hnl[\$)]-1—  a+a  (25) 

and  yields  the  quantum  jump  probability 


l<M0l2 

l<A(0|2 


d*(l-*2) 


=  2 


(26) 


{it+2kYz)  p*(r  5  f  bounds  0 

i  Vi  pi  C r ,  r' bou„d;  t ) 

iV*  [(pi  (r,r;t)  +  (r ,  r ;  t))  pfe  (r,  r'bound;  t )  + 
P%  (r,r;t)  pbJ (r  9  1*  bound  9  01 

iV*  [(p£  (r,r;t)  +  pbb  (r,  r ;  t ))  pi (r ,  r’bound;  t)  + 
pi  (r,  r ;  t)  pi  (r,  r'b0und;  0]  (23) 


Finally,  the  post-jump  wave  function  turns  out  to  be 
e-'T^Vi  4>(t)) 

\<f>(t+m  = - ■  ,  „  (27) 

rwtw 

v  \m2 

i.e.,  <t>e  ( t  +  8 1)  =  0  and 

+ so = e-'?^v-|^(o)/vi^(or/icA(or . 


Here  r'b0und  is  the  position  coordinate  on  the  boundary, 
r'bound  =  (x1,  y',  0)  for  the  forward  propagating  compo¬ 
nent  and  r’bound  =  (x',  y',  L)  for  the  backward 
propagating  component.  Note  that  this  situation  is 
different  from  the  case  of  confined  systems,  where 
Per  ( y,z;y z! bound,  0  =  0  on  the  boundary. 


In  case  of  decay  the  wave  functions  are  given  by 
Eq.  (19). 

These  simulations  are  carried  out  within  the  slowly 
varying  envelope  approximation  introduced  in  Sec.  2. 
The  expressions  for  the  slowly  varying  amplitudes 
Bc(r,r+80.fig(r,r+80,Ef  (rj+ht)  and  Fg(r,r+8f )  can 
be  found  easily  from  the  Eqs.  (17)  and  (19). 


7.  Appendix  B.  Explicit  Expressions 
for  the  Hartree  Monte  Carlo 
Simulations 

The  explicit  form  of  the  incoherent  Liouvillian  Lnc  in 
Eq.  (1)  leads  to  the  system  operators 


c- 


3T 

—  (1  -  AT V-e- 


(24) 


where  3  (1  -  X1)^  is  4ic  dipole  radiation  pattern  for  the 
two-level  atoms  and  x=cos6,  0  being  the  angle  be¬ 
tween  the  quantization  axis  and  the  direction  in  which  the 

spontaneous  photon  is  emitted  and  exp  ^ ±  /  —  e(k )r'j 
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